function [dfTildeOutS_par, domegaTildeIBOutS_par, accPSBT, gyroPSBT, dfCrsDepS, domegaIBCrsDepS, accdKScalCoef, gyrodKScalCoef] ...
    = imuoutdiff_par(IMUParDiff, IMUPar, fB, omegaIBB, omegaIBBRate, iSamp, Tl)
%IMUOUTDIFF_PAR 计算惯组参数误差（微分）引起的惯组输出变化量（微分），即((del vTilde_Out^S)/del p)*dp
%
% Input Arguments
% # fB, omegaIBB, omegaIBBRate: 当Tl不为空时，应传入增量值，此时omegaIBBRate为3*1向量（iSamp-1周期值）
%
% References
% # 《应用导航技术》 误差参数残差补偿算法

if nargin < 7
    Tl = [];
end
if nargin < 6
    iSamp = 1;
end

accPSBT = (IMUPar.acc.PS0B(:, :, iSamp) * (eye(3)+IMUPar.acc.dPSS0(:, :, iSamp)))';
gyroPSBT = (IMUPar.gyro.PS0B(:, :, iSamp) * (eye(3)+IMUPar.gyro.dPSS0(:, :, iSamp)))';

if isempty(Tl)
    dfSizeS = acclevarmerr_cont(accPSBT, IMUPar.acc.lB(:, :, iSamp), omegaIBB(:, iSamp), omegaIBBRate(:, iSamp));
    dfAnisoS = accanisoerr_cont(IMUPar.acc.KAniso(:, iSamp), IMUPar.acc.CBA, omegaIBB(:, iSamp));
    timeMult = 1;
else
    dfSizeS = acclevarmerr_inc(accPSBT, IMUPar.acc.lB(:, :, iSamp), omegaIBB(:, iSamp), omegaIBBRate, Tl);
    dfAnisoS = accanisoerr_inc(IMUPar.acc.KAniso(:, iSamp), IMUPar.acc.CBA, omegaIBB(:, iSamp), omegaIBBRate, Tl);
    timeMult = Tl;
end
dfCrsDepS = dfSizeS + dfAnisoS;
domegaIBGSensS = IMUPar.gyro.DGSens(:, :, iSamp) * fB(:, iSamp);
domegaIBCrsDepS = domegaIBGSensS;

dfBiasExS = IMUPar.acc.dfBiasS(:, iSamp)*timeMult + dfCrsDepS + IMUPar.acc.dfQuantS(:, iSamp)*timeMult;
domegaIBBiasExS = IMUPar.gyro.domegaIBBiasS(:, iSamp)*timeMult + domegaIBCrsDepS + IMUPar.gyro.domegaIBQuantS(:, iSamp)*timeMult;

fTildeScalS = accPSBT*fB(:, iSamp) + IMUPar.acc.dfBiasS(:, iSamp)*timeMult + dfCrsDepS;
omegaTildeIBScalS = gyroPSBT*omegaIBB(:, iSamp) + IMUPar.gyro.domegaIBBiasS(:, iSamp)*timeMult + domegaIBCrsDepS;

accdKScalCoef = scalfacterrsign(IMUPar.acc.dKScalP(:, :, iSamp), IMUPar.acc.dKScalN(:, :, iSamp), fTildeScalS);
gyrodKScalCoef = scalfacterrsign(IMUPar.gyro.dKScalP(:, :, iSamp), IMUPar.gyro.dKScalN(:, :, iSamp), omegaTildeIBScalS);
[~, accKScal] = polyscalfact(fTildeScalS, accdKScalCoef, IMUPar.acc.KScal0(:, iSamp), Tl);
[~, gyroKScal] = polyscalfact(omegaTildeIBScalS, gyrodKScalCoef, IMUPar.gyro.KScal0(:, iSamp), Tl);

% 标度因数误差微分引起的惯组输出变化量，考虑高阶标度因数误差系数
accDiffdKScalCoef = scalfacterrsign(IMUParDiff.acc.dKScalP(:, :, iSamp), IMUParDiff.acc.dKScalN(:, :, iSamp), fTildeScalS);
dfTildeOutS_SF = zeros(3, 1);
for i=1:size(accDiffdKScalCoef, 2)
    dfTildeOutS_SF = dfTildeOutS_SF + accDiffdKScalCoef(:, i).*((fTildeScalS/timeMult).^(i-1));
end
dfTildeOutS_SF = IMUPar.acc.KScal0(:, iSamp) .* dfTildeOutS_SF .* (accPSBT*fB(:, iSamp)+dfBiasExS);
gyroDiffdKScalCoef = scalfacterrsign(IMUParDiff.gyro.dKScalP(:, :, iSamp), IMUParDiff.gyro.dKScalN(:, :, iSamp), omegaTildeIBScalS);
domegaTildeIBOutS_SF = zeros(3, 1);
for i=1:size(gyroDiffdKScalCoef, 2)
    domegaTildeIBOutS_SF = domegaTildeIBOutS_SF + gyroDiffdKScalCoef(:, i).*((omegaTildeIBScalS/timeMult).^(i-1));
end
domegaTildeIBOutS_SF = IMUPar.gyro.KScal0(:, iSamp) .* domegaTildeIBOutS_SF .* (gyroPSBT*omegaIBB(:, iSamp)+domegaIBBiasExS);

% 零偏、标度因数误差、失准角等参数的微分引起的惯组输出变化量
dfTildeOutS_par = diag(accKScal)*IMUParDiff.acc.dfBiasS(:, iSamp)*timeMult + dfTildeOutS_SF ...
    + diag(accKScal)*(IMUPar.acc.PS0B(:, :, iSamp)*IMUParDiff.acc.dPSS0(:, :, iSamp))'*fB(:, iSamp);
domegaTildeIBOutS_par = diag(gyroKScal)*IMUParDiff.gyro.domegaIBBiasS(:, iSamp)*timeMult + domegaTildeIBOutS_SF ...
    + diag(gyroKScal)*(IMUPar.gyro.PS0B(:, :, iSamp)*IMUParDiff.gyro.dPSS0(:, :, iSamp))'*omegaIBB(:, iSamp);
end